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The  results  of  this  research  will  be  submitted  to  an  archival  journal  for  possible 
publication.  Accordingly,  most  of  the  following  sections  of  the  report  are  prepared  in 
the  form  of  such  paper. 

The  results  will  also  be  presented  at  IMECE-2005  (International  Mechanical 
Engineering  Congress)  as  a  part  of  the  symposium  on  adaptive  structures  (invited 
paper). 

Details  of  these  submissions  will  be  forwarded  to  ARO  as  soon  as  the  corresponding 
papers  are  accepted  for  publication  and  presentation. 
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Abstract 


The  paper  is  concerned  with  a  physically  nonlinear  piezoelectric  material  behavior  and  its 
applications  to  practical  problems.  A  survey  of  work  dealing  with  the  phenomenon  is 
included  in  the  introduction.  Subsequently,  the  emphasis  is  on  the  analysis  of  vibrations 
of  piezoelectric  rods  where  a  rather  unique  situation  is  observed,  i.e.  the  response  of  a 
nonlinear  system  can  be  modeled  by  linear  equations  of  motion.  The  solutions  are 
obtained  analytically  by  the  Lagrange  equation  and  by  the  Generalized  Galerkin 
procedure.  Applying  either  of  these  methods,  the  study  of  forced  vibrations  of  a 
physically  nonlinear  piezoelectric  rod  subject  to  a  periodic  electric  field  in  the  axial 
direction  is  reduced  to  the  analysis  of  a  system  of  nonhomogeneous  Mathieu-Hill 
equations.  In  the  particular  case  where  the  interaction  between  axial  and  radial  vibrations 
can  be  neglected,  the  closed-fonn  solution  for  the  former  vibrations  is  obtained  in  the 
paper  and  it  is  shown  that  both  the  Lagrange  equation  and  the  Generalized  Galerkin 
procedure  yield  identical  results.  Numerical  examples  presented  in  the  paper  elucidate 
the  significance  of  physically  nonlinear  effects  that  should  not  be  arbitrary  disregarded  in 
design,  without  a  proper  evidence. 

Introduction 

The  well  known  physically  linear  relationships  between  the  tensors  of  stress  and  strain 
and  the  vectors  of  electric  field  and  electric  displacement  for  piezoelectric  materials  are 
applicable  to  a  particular  case  where  nonlinear  effects  are  negligible.  While  geometric 
nonlinearity  can  be  partially  incorporated  into  these  equations  by  an  appropriate  choice  of 
the  strain-displacement  relationships,  neither  this  nonlinearity  nor  physical  nonlinearity 
of  the  material  are  fully  reflected  in  the  linear  version  of  the  constitutive  equations.  A 
general  form  of  constitutive  equations  accounting  for  nonlinear  products  of  the 
components  of  the  strain  vector  (geometrically  nonlinear  effect)  and  physically  nonlinear 
tenns  has  been  derived  and  published  (see  for  example,  Maugin  et  ah,  1992).  However, 
the  complexity  of  these  equations  as  well  as  the  lack  of  experimental  data  on  the 
coefficients  at  the  nonlinear  terms  and  the  difficulty  involved  in  their  evaluation 
prevented  a  wide  acceptance  of  nonlinear  constitutive  equations  in  design  and  practical 
applications.  This  makes  it  important  to  elucidate  a  relative  contribution  of  geometrically 
and  physically  nonlinear  terms  and  to  assess  both  the  necessity  of  their  incorporation  in 
the  analysis  and  their  relative  qualitative  and  quantitative  effects  on  the  solution. 

Early  studies  of  the  effects  of  stress  and  electric  fields  on  the  response  of  piezoelectric 
materials  were  conducted  by  Berlington  and  Krueger  (1959),  Woolett  and  Leblanc 
(1973),  Krueger  (1954,  1967,  1968a,  1968b),  Brown  and  McMahon,  (1962,  1965)  and 
Fritz  (1978).  A  nonlinear  nature  of  the  problem  is  evident  in  these  investigations.  The 
theoretical  formulation  was  also  developed  by  Tiersten  (1971)  who  later  applied  it  to  the 
problems  of  thin  and  membrane  piezoelectric  plates  subjected  to  high  electric  fields 
(Tiersten,  1993  a,  b).  Other  derivations  of  physically  nonlinear  equations  were  published 
by  Nakagawa  et  al.  (1973)  and  Cho  and  Yamanouchi  (1987). 


5 


Notably,  a  physically  nonlinear  behavior  of  piezoelectric  materials  is  also  reflected  in  the 
hysteresis  or  butterfly  loops  in  the  electric  field  -  strain  or  strain-stress  planes.  Such 
loops  and  relevant  nonlinear  phenomena  were  described  and  discussed  by  a  number  of 
authors  (Chen  and  Montgomery,  1980;  Chen  and  Madsen,  1981;  Bassiouny  et  ah,  1988; 
Leigh  and  Zimmerman,  1991;  Ge  and  Jouaneh,  1995,  etc.) 

Beige  and  Schmidt  (1982)  and  Beige  (1983)  included  higher-order  electric  and  elastic 
tenns  in  their  studies  of  longitudinal  vibrations  of  a  plate  with  a  3 1  piezoelectric  effect. 
This  effect  has  also  been  studied  by  von  Wagner  and  Hagedorn  (2002)  for  piezoelectric 
beams.  Other  studies  of  von  Wagner  (2003,  2004)  and  von  Wagner  and  Hagedorn  (2003) 
were  concerned  with  physical  nonlinearity  in  3 3 -piezoelectrics.  The  work  of 
Chattopadhyay  et  al.  (1999)  employed  nonlinear  constitutive  equations  incorporating  a 
cubic  nonlinearity  for  the  electric  field  to  analyze  helicopter  blades  (modeled  by 
composite  box  beams).  The  subsequent  work  accounted  for  transverse  shear 
deformability  of  monocoque  and  sandwich  plates  combined  with  the  nonlinear 
piezoelectric  effect  of  embedded  or  mounted  sensors  and  actuators  (Thornburgh  and 
Chattopadhyay,  2001). 

Among  recent  studies  that  attempt  to  implicitly  account  for  nonlinearities  by  using 
variable  coefficients  in  linear  constitutive  relationships  one  can  mention  the  paper  by 
Sherritt  et  al.  (1996)  where  the  piezoelectric  coefficient  d33  of  lead  zirconate  titanate 

ceramics  was  shown  dependent  on  stress  as  well  as  being  a  function  of  temperature  and 
frequency.  In  particular,  a  step  stress  caused  a  time  dependent  variation  in  this 
coefficient  that  was  attributed  to  a  slow  movement  of  90°  domain  walls  in  the  material. 
Further  studies  of  the  effects  of  various  factors  on  the  coefficients  in  linear  constitutive 
equations  for  piezoelectrics  were  conducted  by  the  same  group  at  the  Royal  Military 
College  of  Canada.  In  particular,  experimental  data  elucidating  these  phenomena  was 
presented  in  the  papers  of  Wiederick  et  al.  (1996)  and  Sherrit  et  al.  (1996,  1997).  In  these 
studies,  it  was  shown  that  both  the  piezoelectric  constants  as  well  as  the  pennittivity  are 
nonlinear  functions  of  the  applied  electric  field.  The  experimental  approach  employed  to 
measure  piezoelectric  constants  was  based  on  the  so-called  optical  lever  used  to  measure 
the  strains  in  the  1  or  3  directions  as  functions  of  the  applied  electric  field  in  the  3- 
direction.  It  was  also  observed  in  these  studies  that  the  piezoelectric  coefficients  increase 
almost  linearly  with  the  stress.  The  nonlinear  contributions  became  much  more 
pronounced  at  the  applied  electric  field  exceeding  500-1000  V/mm. 

The  intrinsic  and  extrinsic  contributions  to  the  piezoelectric  effect  were  discussed  by 
Yang  et  al.  (2000).  While  the  intrinsic  property  is  typical  for  a  single  domain  crystal,  the 
contributions  associated  with  the  presence  of  multiple  crystals  are  usually  extrinsic.  Such 
contributions  are  dominant  in  soft  piezoceramics  where  they  are  associated  with  domain 
switching  under  the  influence  of  high  electric  fields  (Mukherjee  et  al.,  2001).  In  general, 
it  has  been  observed  that  piezoelectric  coefficients  are  nonlinear  functions  of  the  applied 
compressive  stress.  As  follows  from  this  paper  and  from  other  publications  of  the  same 
group  (see  for  example,  Ren  et  al.,  2000),  piezoelectric  coefficients  d3l,d33,dl5  typically 
increase  with  the  electric  field,  particularly  in  soft  piezoceramics,  while  an  increase  in  the 
frequency  results  in  a  very  small  decrease  of  these  coefficients. 
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Wang  and  Cannan  (1995)  showed  in  their  experiments  that  both  the  coefficient  of 
thennal  expansion  as  well  as  the  piezoelectric  coefficient  J31  are  affected  by  a  cryogenic 

temperature.  Moreover,  it  was  shown  that  J31  is  affected  by  the  magnitude  of  strain.  In 

addition,  the  electric  field  was  recently  shown  to  affect  the  elastic  modulus  (Chaplya  and 
Carman,  2002). 

Barrett  (1995)  showed  a  nonlinear  relationship  between  the  strain  in  the  piezoelectric 
(PZT-5H)  actuator  and  the  applied  electric  field.  The  nonlinearity  was  observed,  even 
though  the  electric  field  considered  in  this  paper  was  limited  to  less  than  500  V/mm. 

Bert  and  Birman  (1998)  proved  that  both  the  stress  and  temperature  affect  the  coefficient 
J3]  in  a  one-dimensional  problem.  In  addition,  they  characterized  the  variations  of  the 

coefficient  of  thermal  expansion  with  the  stress  and  electric  field.  For  example,  it  was 
shown  that  CTE  of  PZT-4  increases  by  10%  if  the  stress  reaches  the  static  strength  value, 
and  even  more  remarkable,  CTE  of  PZT-5A  increased  by  15.8%  under  the  electric  field 
equal  to  2,000  V/mm.  This  work  was  further  expanded  by  Bert  and  Binnan  (1999)  to 
two-dimensional  and  three-dimensional  cases. 

Joshi  (1992)  derived  physically  nonlinear  constitutive  equations  for  piezoceramics  by  the 
assumption  that  material  constants  are  independent  of  the  magnitude  of  stress  or  electric 
field.  This  solution  was  obtained  using  the  thermodynamic  Gibbs  potential  and  retaining 
the  second-order  terms  in  the  total  differentials  of  dependent  variables  (strains,  electric 
flux  density,  and  entropy).  Numerical  examples  were  not  presented  in  this  paper  making 
it  impossible  to  estimate  a  relative  contribution  of  various  terms  in  the  physically 
nonlinear  fonnulation.  An  important  contribution  was  related  to  accounting  for  the 
electrostrictive  and  elastostrictive  effects.  The  constitutive  equations  for  a  physically 
nonlinear  material  derived  by  Maugin  et  al.  (1992)  based  on  the  analysis  of  the  volume 
energy,  i.e.  the  energy  density  of  an  electroelastic  solid,  are  similar  to  equations  of  Joshi 
(1992). 

A  nonlinear  relationship  between  the  deflection  of  the  tip  of  a  bimorph  working  in  the  3 1 
mode  and  the  applied  electric  field  was  observed  for  various  piezoelectric  actuators  by 
Wang  et  al.  (1999)  who  attributed  this  nonlinearity  to  an  increase  of  d3X  with  the  applied 

electric  field.  Note  that  the  physical  nonlinearity  was  reported  in  this  paper,  although  the 
electric  field  was  relatively  low  (150  V/mm). 

The  effect  of  physical  nonlinearity  on  shape  control  of  composite  laminated  beams  was 
considered  by  Achuthan  et  al.  (2001).  As  follows  from  this  study,  the  voltage  required  to 
control  the  shape  of  the  beam  is  significantly  reduced  when  physical  nonlinearity  of 
piezoelectric  patches  on  the  beam  surface  is  taken  in  to  account. 

Notably,  the  electrostrictive  response  can  include  higher-order  nonlinear  tenns,  in 
addition  to  the  well-known  quadratic  relationship  between  the  stress  and  the  electric  field 
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(Sherrit  et  al.,  1999).  However,  these  nonlinear  contributions  become  essential  only  in 
the  case  of  very  high  applied  fields. 

Explicit  expressions  for  piezoelectric  coefficients  dn ,  d33  as  functions  of  the  peak-to- 

peak  voltage  amplitude  were  obtained  for  a  class  of  materials  by  Williams  (2004).  These 
relationships  involved  quadratic  power  of  the  peak-to-peak  voltage,  reflecting  physical 
nonlinearity. 

Representative  empirical  equations  including  physical  nonlinearity  are  shown  below.  For 
example,  the  results  communicated  by  Barrett  (2002)  regarding  the  response  of  G-l  195 
piezoceramics  under  the  electric  field  of  up  to  600  V/nnn  include  the  following  nonlinear 
relationship  between  the  microstrain  and  the  electric  field  (31  effect): 

^  =  (0.227 Ez  +  0.0002437s2  )*  10  6  (1) 

Priya  et  al.  (2001)  found  that  the  relationship  between  the  coefficient  <731  and  the  squared 

applied  elastic  strain  is  linear.  In  particular,  for  soft  PZT-5A,  this  relationship  was 
obtained  in  the  form 

d3X  =2.03*10  10 +1.12*  10^f 2  (2) 

where  the  squared  strain  varied  from  zero  to  2.8*  1  O'7  and  the  piezoelectric  coefficient 
was  measured  in  C/N. 

The  present  paper  concentrates  on  the  investigation  of  physically  nonlinear  effects  on  the 
behavior  of  piezoelectric  rods  polarized  in  the  axial  direction.  Such  problems  are 
important  in  piezoelectric  transducers  (Fig.  1)  used  in  underwater  hydrophones,  acoustic 
imaging  and  medical  applications.  Typically,  the  behavior  of  piezoelectric  rods  and  1-3 
piezocomposites  consisting  of  rods  embedded  in  the  matrix  are  studied  using  linear 
constitutive  equations  (Li  and  Sottos,  1995,  1996a,b,  Sigmund  et  al.,  1998).  Recent 
papers  by  Tan  and  Tong  (2001, 2002)  extended  the  study  to  the  physically  nonlinear 
static  formulation.  Physically  nonlinear  dynamic  problems  were  also  considered  by 
Wagner  (2003,  2004)  and  Wagner  and  Hagedom  (2003).  A  related  problem  of  linear 
vibrations  of  thin  piezoceramic  discs  accounting  for  coupled  axial,  tangential  and  radial 
modes  was  analyzed  by  Huang  et  al  (2004).  In  the  present  paper,  a  unified  formulation  is 
presented  enabling  us  to  identify  essential  physically  nonlinear  effects  and  conduct  a 
comprehensive  analysis  of  the  problem.  As  follows  from  this  paper,  physically  nonlinear 
effects  may  significantly  influence  the  dynamic  response  of  a  piezoelectric  rod. 

Analysis 

1 .  Derivation  of  physically  nonlinear  constitutive  relationships  for  an  orthotropic 
cylindrical  piezoelectric  rod  subject  to  an  electric  field  in  the  axial  direction 


8 


Consider  a  cylindrical  rod  shown  in  Fig.  2  subject  to  an  electric  field  in  the  z-direction. 
For  convenience,  the  following  notations  employed  in  the  subsequent  discussion  are 
introduced  at  this  phase.  The  cylindrical  coordinate  axes  z,  r,  6  identified  in  Fig.  1  are 
also  denoted  as  3,  1  and  2,  respectively.  The  problem  being  axisymmetric,  the 
circumferential  displacement  is  equal  to  zero,  while  the  axial  and  radial  displacements  are 
denoted  by  w  and  u  ,  respectively. 

The  rod  subject  to  a  dynamic  excitation  induced  by  an  electric  field  E .  =  E3  experiences 

vibrations  in  the  axial  and  radial  directions,  while  tangential  motion  occurs  only  if 
axisymmetry  is  violated  by  adjacent  structures.  If  the  side  surface  of  the  rod  is  not 
attached  to  other  structural  elements  of  if  the  adjacent  structure  is  quasi-isotropic,  the  rod 
retains  a  cylindrical  shape  during  these  vibrations.  Therefore,  shearing  stresses  and 
strains  are  equal  to  zero.  In  the  present  analysis,  the  rod  is  assumed  ‘anchored”  at  the 
plane  z  =  0 ,  preventing  axial  displacements  at  this  location. 

The  constitutive  relations  are  derived  following  the  approach  by  Maugin  et  al.  (1992). 
Note  that  notations  used  in  this  book  differ  from  a  number  of  references  where  similar 
relationships  were  derived,  such  as  Joshi  (1992).  However,  accounting  for  this  difference 
does  not  result  in  a  physically  different  formulation. 

The  energy  density  obtained  neglecting  geometric  nonlinearities  is  (Maugin  et  al.,  1992): 


'F  =  —  C„R£  s  n  -  e  E  s  —  —  e  „RE  ££„-  —  £  EE  —  —  s  E  E  E  —  —  l  E  E  £n 

r\  Ctp  (X  fj  HICX  171  (X  r\  WlOCfj  YiX  Ct  fj  YtlYl  HI  tl  WlTip  HZ  VI  p  *  ™  ™  ^ 


mna  m  n  a 


(3) 


where 

Cap  are  elastic  stiffness  constants; 

£a  are  strains; 

ema  are  piezoelectric  constants; 

Em  are  components  of  electric  field  in  the  corresponding  direction; 
emap  are  electroelastic  constants; 

£mn  and  £  are  dielectric  coefficients  (permittivity  and  third-order  dielectric 

coefficients,  respectively); 
l m„a  arc  electrostrictive  coefficients. 


The  constitutive  relations  yielding  stresses  and  electric  displacements  are  obtained  as 


(7 ,,  = 


d£„ 


D„,  =  - 


a^F 

oE,„ 


(4) 
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This  yields  the  following  expressions  for  the  components  of  the  stress  tensor  and  electric 
displacements: 


cr  =  C„a£ a  —  e„„E_  -  „ - E.„EE, 


'aP^P  ^  ma  m  '~mapJ^mc'P 


mna  m  n 


(5) 


D,n  ~  ema£a  +  ~emap£a£p  +  £  nmE  n  +  ~£nmpEnEp  +  lnnaE n£ a 


(6) 


The  terms  underlined  in  the  right  side  of  these  equations  are  retained  in  the  linear 
formulation.  As  indicated  above,  numerous  authors  implicitly  account  for  the 
nonlinearity  by  using  variable  coefficients  of  these  terms,  dependent  on  various  factors. 
Note  that  the  quadratic  function  of  strain  in  (6)  can  usually  be  neglected  even  if  the 
problem  is  geometrically  nonlinear. 


In  the  problem  considered  in  this  paper,  we  are  concerned  with  a  dynamic  response  of  the 
rod  to  the  applied  field  in  the  axial  direction,  i.e.  E3(t)  =  E,{t) .  Accordingly,  the  nonzero 
elements  of  the  tensor  of  stress  given  by  (5)  are 


En  —  e3UE3  EI2  —  e312E3  EI3  ~  e313E3 

>  = 

E12  e32lE3  E22  ~  e322E 3  E 23  ~  e323E3 

< 

CTz. 

Cj3  —  e22lE3  C23  —  e332E^  C33  —  e333fs3^ 

'31 


'32 


if  33 


2 


l331 


‘■332 

1^333 


(7) 


2.  Analysis  of  uncoupled  axial  vibrations 

This  section  presents  the  solution  of  the  problem  of  uncoupled  axial  vibrations  of  a 
piezoelectric  rod  activated  by  a  periodic  in  time  electric  field  in  the  axial  (z)  direction. 
The  solution  is  obtained  by  assumption  that  the  effect  of  radial  vibrations  on  the  axial 
motion  can  be  disregarded. 

Two  methods  of  the  analysis  are  considered,  namely  the  Generalized  Galerkin  procedure 
and  the  Lagrange  equation.  In  particular,  the  fonner  procedure  is  different  from  the 
conventional  and  well  known  Galerkin  procedure  in  enabling  one  to  address  the  issue  of 
boundary  conditions  that  cannot  be  satisfied  by  expressions  chosen  to  represent 
displacements.  Although  this  procedure  is  relatively  little  known,  the  solutions  utilizing 
it  have  been  published  (Simitses,  1986;  Houbolt  and  Brooks,  1958).  The  choice  of  the 
appropriate  sign  for  the  term  incorporating  the  boundary  conditions  into  the  Generalized 
Galerkin  fonnulation  can  only  be  established  if  the  equations  of  equilibrium  and  the 
boundary  conditions  are  derived  from  the  Hamilton  principle.  Accordingly,  the  analysis 
begins  with  the  consideration  of  this  principle  in  application  to  the  present  problem. 

2.1.  Derivation  of  the  equation  employed  in  the  Generalized  Galerkin  procedure  from  the 
Hamilton  principle 
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The  Hamilton  principle  applied  to  the  problem  of  axial  vibrations  of  a  piezoelectric  rod 
where  the  motion  is  excited  by  an  electric  field  and  external  mechanical  forces  are  absent 
is  formulated  as 


I'1  (SV  -  Sr)dt  =  0  (8) 


where  the  integration  is  conducted  over  an  arbitrary  interval  of  time  t0  <t  <tx,  V  is  the 
strain  energy,  and  T  is  the  kinetic  energy. 

The  strain  energy  accumulated  in  the  rod  that  experiences  vibrations  in  the  axial  direction 
is  given  by 


2  n  a  h 

SV  =  111  <7  _  §£_  r dzdrd  6 
0  0  0 


(9) 


where  a  is  the  radius  of  the  rod,  and  the  variation  of  the  axial  strain  is  related  to  the 
variation  of  the  axial  displacement  w  by 


8s :  =  (Sw\: 


The  substitution  of  (10)  into  (9)  and  the  integration  by  parts  yields 


8V  =  m2j  {crz8w\=h  -  ( a:Sw)!=0  -  J (crz ),z  dwdz 
^0  -  0 

The  variation  of  the  kinetic  energy  is 

h 

8T  =  m2mj  wSwdz 
o 

where  m  is  the  mass  density  of  the  rod  material. 
Integrating  this  expression  by  parts  one  obtains 


\dt 


n 

ST  =  m2m | 


ll 

(wSw\=ti  ~(wSw)t=to  -jwSwdt 


dz 


(10) 


(ID 


(12) 


(13) 


The  first  two  terms  under  the  integral  in  the  right  side  of  (13)  can  be  taken  equal  to  zero 
assuming  that  8w(t  =  tl)  =  8w(t  =  /0 )  =  0  (Whitney,  1987).  Then  substituting  the 

variations  of  strain  and  kinetic  energies  given  by  (1 1)  and  (13),  respectively,  into  (8)  one 
obtains  both  the  equation  of  motion  as  well  as  the  boundary  conditions  that  are  actually 
well  known  (the  reason  for  the  previous  derivation  is  explained  below): 
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cr.  =  0  or  Sw  =  0 
a,  =  0  or  Sw  =  0 


(14) 


mw  -  cr.,z  =  0 

z  =  0 : 
z  =  h: 


(15) 


The  boundary  conditions  that  can  be  specified  from  (15)  reflect  that  in  the  present 
problem  the  cross  section  z  =  0  is  prevented  from  axial  motion,  so  that c5W(z  =  0)  =  0 . 
The  axial  stresses  at  the  free  end  of  the  rod  should  be  equal  to  zero,  i.e.  a,(z  =  h)=  0 . 

In  the  case  where  it  is  impossible  to  satisfy  the  stress  boundary  conditions,  following  the 
Generalized  Galerkin  procedure,  the  displacements  can  be  sought  in  the  form 

MzS)  =  XZ’(z)r’(0  (16) 

i 


and  the  system  of  equations  of  motion  with  respect  to  the  functions  of  time  Tj  (l  ) 
available  from  (8)  by  substituting  (1 1)  and  (13)  is 

h 

J (<JZ , z  —mw)Zi (z)dz  -  cr_(z  -  h)Zt(h)  =  0  (17) 

o 


It  is  observed  that  the  sign  of  the  second  term  in  (17)  would  not  be  evident  without  the 
previous  derivation. 

2.2  Analysis  using  normal  modes  to  represent  the  motion 

The  motion  considered  in  this  section  is  represented  in  terms  of  normal  modes  of  axial 
vibrations  of  the  rod.  Two  methods  employed  to  analyze  the  problem  are  the  Generalized 
Galerkin  procedure  and  the  Rayleigh-Ritz  method.  As  is  shown  below,  these  methods 
yield  identical  systems  of  equations  of  motion. 

The  axial  free  vibrations  of  a  rod  with  the  boundary  conditions  specified  above  are 
represented  in  the  form 


w(x,  t)  =  ^  (Au  sin  Aft  +  A2i  cos  /fy)sin 


17TZ 

2 Ti 


(18) 


where  i  is  a  natural  odd  number,  Akt  are  constants  of  integration  specified  from  the 


in 


initial  conditions,  and  A.  =  — 


C. 


33 


2  h  V  m 


is  a  natural  frequency. 


2.2.1.  Generalized  Galerkin  procedure 
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Following  the  Galerkin  procedure,  the  axial  displacements  should  be  chosen  in  the  form 
of  series  satisfying  the  boundary  conditions.  However,  it  appears  impossible  to  choose 
such  series,  while  satisfying  the  conditions  both  at  the  “anchored”  cross  section  z  =  0  as 
well  as  at  the  free  end  of  the  rodz  =  h .  Therefore,  the  solution  is  sought  in  the  form 

w=S^(?)/2sin^f  09) 


The  substitution  of  (19)  into  (17)  where  Z(.  =  sin —  yields  a  system  of  uncoupled 

2  h 


equations.  In  particular,  the  i-th  equation  of  this  system  is 


mh 1 


,•2  2 
i  n 


^+(C„-emEz)—Wi  = 


E2  +' 


-/  E 

1 333^2 


sin- 


in 


(20) 


2.2.2.  Solution  by  the  Lagrange  equation 

The  energy  density  within  a  rod  vibrating  in  the  axial  direction  can  be  obtained  from  (3): 


1  2  j-,  1  2I  2I  3I  2 

2  33  z.  33  z.  z  r\  333  z  z  33  z  333  z  333  z  z 

2  2  6  2 


(21) 


Note  that  the  terms  dependent  only  on  the  electric  field  that  are  underlined  in  (21)  do  not 
affect  the  subsequent  solution  and  are  omitted. 

If  axial  vibrations  of  the  rod  are  represented  by  series  (19),  the  Lagrange  equation  for  the 
n-th  tenn  of  these  series  is 


d_ 

dt 


r  dT_  A 


dT  dV 
dW„  +  dW„ 


where  the  kinetic  and  strain  energy  contributions  are  given  by 


T  = 


2 


2 7ih  a 

J  J  J  mw2rdrdzdO 

000 


and 


2 7ih  a 

y=\\\  WrdrdzdO 
000 


(22) 


(23) 


(24) 


respectively. 
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The  substitution  of  (19)  into  (21,  23,  and  24)  and  the  subsequent  application  of  the 
Lagrange  equation  (22)  yields  the  same  system  of  equations  of  motion  as  that  obtained  by 
the  Generalized  Galerkin  procedure,  i.e.  (20). 

2.3  Analysis  using  power  series  to  represent  the  motion 

The  solution  can  also  be  sought  in  the  power  series,  i.e. 

w=£fVj(ty  (25) 

7=1 

where  W j  (/)  are  unknown  functions  of  time  ( j  are  natural  numbers  that  can  be  both  even 
and  odd). 

2.3.1.  Generalized  Galerkin  procedure 

It  is  immediately  obvious  that  series  (25)  satisfy  the  boundary  condition  w{z  =  0)  =  0 , 
while  the  condition  of  zero  axial  stress  at  the  free  end  is  violated.  Accordingly,  the 
Generalized  Galerkin  procedure  implies  that  the  n-th  equation  of  the  system  of  equations 
of  motion  is 


h 

J  (cr„ -mw)z"dz  -  cj,  (z  =  h  )h "  =  0 

o 


(26) 


Note  that  (26)  was  written  accounting  for  the  axisymmetry  of  the  problem  and 
accordingly,  avoiding  the  integration  in  the  circumferential  direction. 


The  substitution  of  the  axial  stress  from  (7)  and  axial  displacements  presented  by  (25) 
into  (26)  and  using  the  strain-displacement  relationship  s,  =  w,_  yields  the  system  of 
coupled  equations  of  motion.  The  n-th  equation  is 


m 


I 


h 


=i  J  +  n  + 


7+«+i  —  f  _  ih^+n~l  _ 

~W,  +(C„ -e,„Ei  h"Wl  +  »Z7 - ~Wj 

t  l  j=2j  +  n-\ 


eee  E-  +  ^133  ^ z 

2  J 


h"  (27) 


2.3.2.  Solution  by  the  Lagrange  equation 

The  energy  density  within  a  rod  vibrating  in  the  axial  direction  is  given  by  (21).  If  axial 
vibrations  of  the  rod  are  represented  by  the  power  series  (25),  the  Lagrange  equation 
yields  the  following  system  of  equations 


j+n+l  _ 


m 


hj+n~i  — 


e33 


V 


}_ 

2  '  333 J 


h 

J 
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(28) 

It  is  easy  to  show  that  the  systems  of  equations  (27)  and  (28)  are  identical,  i.e.  the 
Generalized  Galerkin  procedure  and  the  Lagrange  equation  yield  identical  equations  of 
motion,  similar  to  the  case  where  the  mode  shape  of  vibration  was  represented  by  normal 
modes. 


2.4,  Solution  of  equations  of  axial  motion 

The  systems  of  equations  (20)  and  (28)  can  be  integrated  if  the  electric  field  is  a  known 
function  of  time.  In  the  present  work  it  is  assumed  that  this  function  is  periodic,  i.e. 

E:(t)  =  E  cos  cat  (29) 


Then  the  substitution  of  (29)  into  one  of  the  above-mentioned  systems  of  equations  yields 
a  system  of  equations  of  motion  that  can  be  solved  either  analytically  or  numerically.  For 
example,  the  system  (28)  becomes 


h 


j+n+ 1  _ 


mV - W :  +  (C„  -  e„,E  cos  cot)nY'  i 

TU  +  n  + 1  '  V  33  333  ’  i 


hi+n~l  — 


j=\  J  +  n  + 

'\  ,  1 

—  lmE  +  e^E  cos  cot  +  —  lmE ~  cos  loot 


~t  j  +  n- 1 


(30) 


hn 


The  substitution  of  (29)  into  (20)  yields 


■  2  2 

j  n 


^-W,  +(C„-emEcosa,t)^-Wl  = 


1  2  1  2 

— l333E  +  e33E  cos  cot  +  —  lmE  cos2  cot 


sin- 


in 


(31) 


The  systems  of  equations  (30)  and  (31)  represent  systems  of  nonhomogeneous  Mathieu- 
Hill  equations.  Such  equations  were  encountered  since  the  first  part  of  the  last  century 
(Strutt,  1932;  McLachlan,  1947)  when  studies  were  often  limited  to  homogenous 
equations  that  are  employed  to  solve  problems  of  dynamic  or  parametric  stability  of 
structures  (Bolotin,  1964).  A  typical  problem  where  the  analysis  is  reduced  to  a 
nonhomogeneous  Mathieu  equation  is  the  motion  of  a  rod  with  initial  imperfection 
subject  to  a  periodic  in  time  axial  force  (Bolotin,  1964).  A  related  problem  is  that  of 
forced  vibrations  of  a  rod  subject  to  an  eccentrically  applied  driving  force.  In  addition  to 
the  classical  monograph  of  Bolotin  (1964),  a  number  of  investigations  have  been 
concerned  with  the  problem  of  an  interaction  between  forced  and  parametric  vibrations 
(Hsu  and  Cheng,  1974;  Nguyen,  1975;  Troger  and  Hsu,  1977;  HaQuang  et  al.,  1987a, b; 
Plaut  et  al.,  1990;  Nguyen  and  Ginsberg,  2001). 


The  analytical  solution  of  the  system  of  equations  (30)  can  be  sought  in  the  form  of 
trigonometric  time  series: 
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(32) 


T-j 7  .  V1  .  run  v'  n 

W ,  =  A.  +  >  A  .  cos - h  >  5  ,  sin  — — 

J  07  ^  rJ  ry  ^  PJ  ry 

r  ^  p  ^ 

where  A.,  A-  and  Bpj  are  unknown  coefficients.  These  coefficients  can  be  determined 

by  substituting  the  series  (32)  into  (30)  or  (31),  equating  the  coefficients  at  the  same 
trigonometric  functions  of  time  and  solving  the  system  of  resulting  linear  algebraic 
equations  with  respect  to  A0j ,  A.  and  B  pj .  The  solution  of  (3 1 )  can  be  obtained  by  the 

same  approach.  Details  of  the  straightforward  solution  of  the  system  of  linear  equations 
with  respect  to  A0j ,  Arj  and  B pj  are  omitted  here  for  brevity. 


If  the  axial  displacement  is  assumed  to  be  a  linear  function  of  the  distance  from  the  cross 
section  z  =  0 ,  i.e.  n  =  j  =  1 ,  and  the  electric  field  is  given  by  (29),  the  system  (30)  is 
reduced  to  a  single  equation 


mh  2  . .  _  _  1^  1  , 

- Wx  +  C33FFt  -e333EWl  cos  cot  =  —l333E2  +  e33 E  cos  cot  +  —l333E2  cos2cot  ^3) 


The  corresponding  equation  for  the  solution  representing  the  motion  in  terms  of  normal 
modes  is  available  from  (31): 


^2  7[2  \  1 

Wx  +(C33  -e333E  coscot)— Wx  =—I333E 2  +  e33E  cos  cot  +  —l 333E 2  cos2  cot  (34) 


As  is  shown  below,  single-term  solutions  that  result  in  equations  (33)  and  (34)  are  always 
accurate  if  the  driving  frequency  co  is  smaller  than  the  fundamental  frequency  of  the  rod. 

The  systems  of  equations  (30)  and  (31)  or  equations  (33)  and  (34)  could  be  simplified  in 
case  of  numerous  piezoelectric  materials  as  shown  below.  It  can  be  observed  that 
physical  nonlinearity  in  the  expression  for  the  energy  density  as  given  by  (21)  results  in 
terms  proportional  to  e333  and  l333  .  Consider  for  example,  equation  (34)  where  it  is 
instructive  to  compare  the  magnitude  of  the  following  coefficients: 

k1  =  e333EW1 

k2=EmE2  (35) 

k3  =  e33E 


Note  that  a  physically  linear  formulation  can  be  obtained  by  setting  k{  =  k2  =  0  . 

Two  materials  chosen  for  the  following  comparison  and  considered  in  numerical 
examples  are  PZT-5H  and  PZN-4.5%PT  (Tan  and  Tong,  2001,  2002).  In  addition, 
LiNbCb  (Maugin,  et  ah,  1992)  is  used  in  the  analysis  of  the  coefficients  (35).  The 
magnitude  of  typical  electric  fields  usually  varies  from  zero  to  2.0  MV/m  (Maugin  et  ah, 
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1992).  The  values  of  material  constants  needed  to  evaluate  the  coefficients  in  (35)  are 
listed  in  Table  1. 

Table  1.  Material  constant  of  representative  piezoelectric  materials 


Material 

*»  (fj 

e333  l42 ) 

C 33  ( GPa ) 

/  (F/  ) 

1 333  \/m> 

PZT-5H 

24.54  (or  3 1.02) 

5.7*1 04 

108.0 

-21.21*1 0"6 

PZN-4.5%PT 

12.1  (or  13.02) 

4.2*  104 

89.0 

- 1 .6 1  *  10'6 

LiNb03 

1.3 

-17.3 

24.5 

-2.76*  10'9 

Using  a  high  electric  field,  i.e.  2.0  MV/m,  since  it  results  in  a  stronger  nonlinear  effect, 
one  obtains  the  values  of  the  coefficients  in  (35)  listed  in  Table  2. 

Table  2.  Coefficients  of  the  equations  of  motion  for  representative  materials 


Material 

K 

k2 

k2 

PZT-5H 

11.4*10luWi 

-2.12*1 07 

5.09*107 

PZN-4.5%PT 

8.9*10luWi 

- 1 .6 1  *  106 

2.42*  107 

LiNbCb 

-34.6*106Wi 

-2.76*103 

2.6*1 06 

3 

Note  that  Wt  in  Table  1  is  nondimensional  and  the  units  of  all  coefficients  are  C*V/m  . 
A  comparison  of  the  coefficients  in  Table  2  yields  the  conclusion  that  the  terms 
proportional  to  the  squared  electric  field,  i.e.  k2 ,  can  be  neglected  if  the  electric  field 
remains  within  certain  limits.  In  particular,  these  terms  can  be  neglected  for  PZN- 
4.5%PT  and  for  LiNb03  even  at  2.0MV/m,  while  the  term  proportional  to  k2  can  be 
neglected  for  PZT-5H  if  the  field  remains  smaller  than  0.5MV/m  since  at  this  electric 
field  k2  is  an  order  of  magnitude  higher  than  k: .  Accordingly,  if  the  electric  field  is 

below  the  limits  specified  above,  the  systems  of  equation  (30)  and  (3 1)  as  well  as 
equations  (33)  and  (34)  can  be  simplified.  In  particular,  the  latter  equations  are  reduced 
to  a  nonhomogeneous  Mathieu  equation: 

Wx  +  22W1  -  2 ju21Wl  cos  cot  =  p  cos  cot  (36) 


where  the  fundamental  frequency,  the  parametric  loading  coefficient  and  the  forcing 
function  are 


1 

3  r 

e333E 

3e33E 

~Ji\ 

m 

~  2 C33 

Pp 

mh 2 

n 

c33 

Mn 

e333E 

n 

2e33E 

“2 Ti 

V  m 

~  2 C33 

P  n 

mh 2 

for  (33) 

(37) 

for  (34) 


where  the  subscripts  identify  the  solutions  in  power  series  (p)  and  in  nonnal  modes  (n). 
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As  will  be  shown  below,  viscous  damping  is  essential  in  the  case  where  the  driving 
frequency  is  close  to  the  fundamental  frequency  of  the  rod.  The  corresponding  expansion 
of  (36)  is 

Wx  +  J3WX  +  A2WX  -2juX 2WX  cos  cot  =  p  cos  cot  (38) 

where  (3  =  ,  c  being  a  damping  coefficient  (for  both  eqns.  33  and  34).  The  value  of 

the  damping  coefficient  can  be  evaluated  from  published  data  for  a  quality  factor  Q  as 

c  =  2  ccrQ  (39) 

where  ccr  is  a  critical  damping  coefficient  of  the  rod  experiencing  axial  vibrations. 

Equations  (36)  or  (38)  are  similar  to  the  equation  for  lateral  vibrations  of  an  imperfect  rod 
subjected  to  a  periodic  in  time  axial  force  (Mettler,  1941;  Bolotin,  1964).  Following 
these  references,  the  solution  for  the  steady  state  vibrations  in  the  vicinity  of  the 
secondary  region  of  parametric  instability,  i.e.  in  the  case  where  the  frequency  of  the 
electric  field  is  close  to  the  fundamental  frequency  can  be  adequately  predicted  retaining 
only  three  terms  in  series  (32): 

Wx  =  A0  +  A2  cos  cot  +  B2  sin  cot  (40) 

Using  (40),  the  ratio  of  the  amplitude  of  vibrations  neglecting  physical  nonlinearity,  i.e., 
using  /j  =  0  in  (36),  to  the  corresponding  amplitude  accounting  for  the  physically 
nonlinear  effect  is  obtained  in  the  fonn 


In  the  presence  of  viscous  damping,  this  ratio  is 


3.  Solution  for  coupled  axial-radial  axisymmetric  vibrations  by  the  Generalized 
Galerkin  procedure 

Equations  of  the  axisymmetric  motion  of  the  rod  are 


18 


=  mu 


(43) 


do,.  |  ar-ae 
dr  r 


do, 

— -  =  mw 
dz 


where  the  terms  in  the  right  side  represent  inertias  in  the  radial  and  axial  directions. 
The  solution  must  satisfy  the  following  conditions: 


r  =  0  :  u  =  0 

z  =  0 :  w  =  0 

r  =  a:  cr, .  =  0 

z  =  h:  crz  =  0 


(44) 


It  is  evident  that  the  exact  satisfaction  of  the  static  (stress)  conditions  in  (44)  is 
impossible. 

The  following  approach  to  the  solution  utilizes  the  Generalized  Galerkin  procedure.  For 
example,  the  displacements  can  be  sought  in  the  form  of  power  series  that  satisfy  the  first 
two  conditions  (44): 

u  =  YJUJ(t)ri  w  =  Y,Wj{t)zj  (45) 

j= i  j= i 


The  Generalized  Galerkin  procedure  yields  a  set  of  equations.  In  particular,  the  n-th 
equation  of  motion  in  the  radial  direction  (j=n)  is 


2 7th  a 

fill 

0  0  0 


for 

dr 


-  +  - 


cr  ~  <Ja 


-  -  mu 


2  nh 


*"+1  dr dzd 6  -  J  Jcr .  (r  =  a)a'l+ldzdd  =  0 


o  o 


(46) 


The  n-th  equation  of  the  axial  motion  obtained  by  this  method  becomes 


2 nh  a 

ill 


0  0  0 


mw 


rz 


2  7T  a 

”  dr  dzd  9  -  J  Jcr.  (z  =  h)rh"  dr d  9  =  0 
o  o 


(47) 

The  substitution  of  the  series  (45)  into  the  constitutive  relations  (7)  and  the  subsequent 
use  of  the  Generalized  Galerkin  procedure  (46)  and  (47)  yield  a  system  of  coupled  time- 
dependent  differential  equations  for  U -(t)  and  IT  (?) . 


In  the  case  where  the  solution  is  sought  using  nonnal  modes  of  motion  as  generalized 
coordinates,  series  (45)  are  replaced  with 
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(48) 


=  Zc/.(Osin 


1=1 


i  Tu¬ 
la 


w 


=  Z^(')s  in 


i'=l 
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where  i  is  an  odd  number.  In  this  case,  the  Generalized  Galerkin  procedure  implies 
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If  the  driving  electric  field  is  a  periodic  and  harmonic  function  of  time  as  is  assumed  in 
this  paper  (see  Eqn.  29),  the  solution  of  the  system  of  time-dependent  equations  of 
motion  available  from  (46),  (47)  or  from  (48),  (49)  can  be  sought  in  the  fonn  of  series 
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Numerical  results  and  discussion 


The  following  results  were  obtained  for  two  representative  materials,  i.e.  PZT-5H  and 
PZN-4.5%PT,  using  data  from  Table  1  (the  values  of  e33  shown  without  brackets  were 

employed  in  calculations).  In  the  examples  presented  below  the  amplitude  of  the  electric 
field  was  chosen  equal  to  2.0MV/m  and  0.5MV/m  for  PZN-4.5%PT  and  PZT-5H, 
respectively.  Accordingly,  the  tenns  proportional  to  the  coefficient  k2  could  be  neglected. 
Notably,  even  in  the  case  of  a  relatively  low  electric  field  for  PZT-5H  the  effect  of 
physical  nonlinearity  was  significant  as  is  shown  in  the  following  examples. 

The  effect  of  physical  nonlinearity  on  the  accuracy  of  the  prediction  of  vibrations  of  a 
piezoelectric  rod  subject  to  electric  fields  with  the  driving  frequency  close  to  the 
fundamental  frequency  of  the  rod  is  illustrated  in  Figs.  3  and  4.  The  horizontal  axes  in 

these  figures  represents  =  ,  so  that  in  the  absence  of  damping  the  amplitude  ratio 

becomes  infinite  at  the  resonant  frequency  according  to  (41).  As  follows  from  Figs.  3 
and  4,  neglecting  physical  nonlinearity  may  cause  a  significant  numerical  error  in  the 
vicinity  of  the  fundamental  frequency,  particularly  at  high  electric  fields.  However,  even 
at  a  relatively  low  electric  field  (Fig.  3),  quantitative  differences  between  physically 
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nonlinear  and  linear  solutions  cannot  be  disregarded  if  F  — >  1 .  It  is  noted  that  each  of  the 
curves  in  Figs.  3  and  4  has  extreme  values  R  =  0  and  R  =  <x> ,  but  some  of  these  values 
could  not  be  reflected  in  the  figures  due  to  a  difficulty  in  graphing  the  results  for  the 
corresponding  pinpointed  values  of  F  .  Nevertheless,  all  essential  trends  are  accurately 
identified  in  these  figures. 

The  effect  of  damping  on  the  ratio  of  amplitudes  of  linear  and  nonlinear  vibrations  is 
elucidated  in  Figs.  5  and  6.  The  quality  factor  employed  to  generate  these  results  was 
chosen  based  on  published  data  for  PZT-5H.  This  factor  was  also  applied  to  PZN- 
4.5%PT  since  damping  data  for  this  material  is  not  available.  However,  independent  of 
the  exact  damping  coefficient  value,  it  is  evident  that  the  qualitative  effect  of  damping  is 
present  only  in  the  immediate  vicinity  of  the  fundamental  frequency  (the  same  conclusion 
is  valid  for  other  resonant  frequencies,  though  the  effect  becomes  weaker  at  higher  modes 
of  motion).  The  effect  of  damping  on  vibrations  of  physically  nonlinear  piezoelectric 
rods  is  further  elucidated  in  Figs.  7  and  8.  As  follows  from  these  figures,  the  amplitudes 
of  vibration  are  almost  unaffected  by  damping  outside  a  very  narrow  spectrum  of  driving 
frequencies  encompassing  the  fundamental  frequency.  However,  at  the  resonant 
frequency  damping  reduces  the  ratio  of  the  amplitude  of  nonlinear  vibrations  with 
damping  to  that  without  damping  (denoted  RD  in  Figs.  7  and  8)  to  zero.  This  is 
anticipated  since  the  amplitude  of  undamped  motion  is  infinite  at  the  resonance. 

The  previous  results  were  generated  using  the  nonnal  mode  approach.  Therefore,  it  is 
important  to  compare  the  solutions  obtained  by  this  method  to  the  results  generated  by 
the  power  series  approach.  Such  comparison  is  shown  in  Figs.  9  and  10  where  the  ratio 
Rpq  represents  a  ratio  of  the  amplitudes  of  nonlinear  vibrations  obtained  by  the  power 
series  and  nonnal  modes  approaches.  The  interpretation  of  the  results  for  nonlinear 
vibrations  shown  in  these  figures  should  account  for  the  fact  that  two  driving  frequencies, 
one  of  them  yielding  a  ratio  Rpn=0  and  the  second  corresponding  to  an  infinite  Rpn 
reflect  the  resonances  of  the  rod  by  the  corresponding  solutions.  Naturally,  the  ratios  Rpn 
exhibit  abrupt  variations  within  the  resonance  region.  Outside  the  resonance  region,  the 
solutions  are  reasonably  close  to  each  other  but  the  power  series  approach  requires 
retention  of  more  tenns  to  accurately  represent  the  motion. 

The  limits  of  the  accuracy  of  a  one -term  solution  by  the  normal  mode  method  are 
elucidated  in  Figs.  1 1  and  12.  As  follows  from  these  figures,  the  contribution  of  higher 
modes  is  negligible  if  the  driving  frequency  is  smaller  than  the  fundamental  frequency  of 
the  rod.  As  the  driving  frequency  increases,  the  contribution  of  the  second  normal  mode 
becomes  essential  (it  was  shown  that  higher  modes  have  a  negligible  effect  in  the  range 
of  frequencies  considered  in  Figs.  1 1  and  12).  The  increase  of  the  contribution  of  the 
second  mode  is  particularly  pronounced  in  case  of  a  high  electric  field. 

Finally,  the  multi-mode  solutions  for  the  nondimensional  amplitude  of  physically 
nonlinear  axial  vibrations  of  piezoelectric  rods  within  a  broad  range  of  driving 
frequencies  are  illustrated  in  Figs.  13,  14  and  15.  As  is  shown  in  these  figures,  the 
amplitudes  of  motion  remain  quite  small,  except  for  very  narrow  regions  encompassing 
the  resonant  frequencies.  However,  though  geometric  nonlinearity  can  be  neglected 
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(with  a  possible  exception  of  the  above-mentioned  regions),  the  physical  nonlinearity 
should  be  accounted  for  as  follows  from  the  previous  discussion. 

Conclusions 

The  effect  of  physical  nonlinearity  on  the  dynamic  response  of  piezoelectric  rods 
polarized  in  the  axial  direction  and  subject  to  a  hannonic  electric  field  acting  in  this 
direction  has  been  investigated.  Such  problems  are  important  in  a  number  of  applications 
employing  piezoelectric  transducers.  It  is  shown  that  physically  nonlinear  effects  become 
more  pronounced  at  higher  fields.  Mathematically,  the  presence  of  physical  nonlinearity 
results  in  additional  terms  in  the  system  of  equations  of  motion,  so  that  these  equations 
change  from  the  equations  for  forced  vibrations  in  the  physically  linear  case  to  the  mixed 
forced-parametric  dynamic  formulation  (Mathieu-Hill  equations). 

Physically  nonlinear  vibrations  considered  in  the  paper  were  investigated  using  the 
Generalized  Galerkin  procedure  and  the  Lagrange  equation.  The  former  formulation 
enabling  us  to  incorporate  all  boundary  conditions  was  derived  from  the  Hamilton 
principle.  The  solution  was  developed  using  two  different  systems  of  generalized 
coordinates,  presenting  the  motion  in  terms  of  nonnal  modes  of  free  vibrations  and  in 
power  series  of  the  coordinates.  Remarkably,  in  each  of  these  cases,  the  Generalized 
Galerkin  procedure  and  the  Lagrange  equation  yield  identical  equations  of  motion.  The 
solution  considered  in  this  paper  was  confined  to  uncoupled  axial  vibrations.  However, 
the  approach  to  the  solution  of  a  coupled  axial-radial  axisymmertic  vibration  problem  has 
also  been  outlined. 

Vibrations  of  the  rods  driven  by  a  harmonic  electric  field  were  numerically  investigated 
in  the  vicinity  of  the  fundamental  frequency.  As  follows  from  the  representative 
examples,  the  error  due  to  neglecting  physical  nonlinearity  may  be  significant.  This  error 
increases  in  the  close  vicinity  to  the  resonant  frequency.  The  inaccuracy  of  the  linear 
analysis  becomes  larger  at  higher  electric  fields. 

Damping  has  a  noticeable  effect  on  vibrations  of  piezoelectric  rods  in  the  vicinity  of  the 
resonant  frequency.  However,  the  frequency  range  where  the  effect  of  damping  is 
significant  is  very  narrow.  Within  this  range,  damping  results  in  a  finite  amplitude  of 
vibrations,  while  without  damping,  the  amplitude  is  infinite,  similarly  to  the  situation 
encountered  in  the  problem  of  forced  vibrations. 

The  analysis  of  vibrations  based  on  using  normal  modes  as  generalized  coordinates 
illustrates  that  higher  modes  have  a  negligible  effect  at  the  driving  frequencies  that  are 
smaller  than  or  close  to  the  fundamental  frequency  of  the  rod.  If  the  driving  frequency 
increases  beyond  the  fundamental  frequency,  the  effect  of  higher  modes  becomes  more 
pronounced,  particularly  at  high  electric  fields.  Therefore,  a  single-mode  analysis  can  be 
safely  employed  only  if  the  driving  frequency  remains  below  or  close  to  the  fundamental 
frequency,  while  a  multi-mode  analysis  is  recommended  at  higher  values  of  the  fonner 
frequency.  The  multi-mode  analysis  conducted  for  axial  vibrations  in  the  range  of 
driving  frequencies  including  the  first  two  natural  frequencies  illustrated  the  presence  of 
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clearly  identified  peaks  in  the  vicinity  of  the  fundamental  and  second  resonant 
frequencies.  Outside  these  peaks,  the  amplitude  of  motion  remained  small,  even  at  high 
electric  fields. 

A  comparison  between  the  results  generated  using  normal  modes  and  the  solution 
employing  power  series  to  model  the  motion  illustrated  that  these  two  solutions  predict  a 
qualitatively  similar  behavior.  The  exception  was  found  in  the  vicinity  of  resonant 
frequencies  where  the  solutions,  while  qualitatively  similar,  may  predict  significantly 
different  amplitudes  at  the  same  driving  frequency.  In  general,  it  is  anticipated  that  the 
nonnal  mode  solution  is  more  accurate,  but  as  the  number  of  tenns  in  the  series 
representing  the  motion  increases,  the  difference  between  the  solutions  should  become 
smaller. 

Note  that  the  effect  of physical  nonlinearity  was  also  investigated  as  a  part  of  the 
program  on  the  development  of  Thunder  actuators  for  QorTek,  Inc.  for  control  surfaces 
of  minituarized  munitions  conducted  by  Dr.  Birman.  The  position  of  the  control  surface 
was  changed  by  activating  piezoelectric  layers  that  caused  an  appropriate  bending  of  this 
surface.  As  was  shown  analytically  (and  confirmed  in  experiments),  accounting  for 
physical  nonlinearity  of  piezoelectric  actuator-layers  resulted  in  a  significant  different 
response  compared  to  a  simplified  approach  neglecting  such  phenomenon.  The 
comparison  between  physically  nonlinear  and  experimental  results  for  representative 
control  surfaces  was  favorable  (see  details  in  Appendix). 
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Appendix:  Physically  nonlinear  effects  in  Thunder  bimorphs  used  in  control  surfaces  of 
miniature  munitions  (development  work  of  Dr.  Victor  Birman  conducted  for  QorTek, 

Inc). 

The  following  brief  description  refers  to  a  program  on  the  development  of  miniaturized 
munitions  where  Dr.  Binnan  conducted  the  design  and  development  work  on  bimorph 
Thunder  actuators.  While  details  of  the  program  cannot  be  disclosed  due  to  the 
proprietary  agreement,  some  results  are  shown  here. 

The  Thunder  actuators  considered  in  the  program  represented  clamped  surfaces  with  a 
central  substrate  bounded  by  piezoelectric  layers  that  were  in  turn  covered  by  thin  cover 
sheets  (Pig.  16).  The  electric  field  applied  in  the  z-direction  resulted  in  the  appropriate 
strains  in  the  longitudinal  direction  (31  effect).  The  results  for  on  of  two  designs 
considered  in  the  project  are  shown  in  Pig.  17.  A  different  length  (L)  considered  in  the 
computations  reflected  uncertainty  about  the  unsupported  length  of  the  cantilevered  beam 
clamped  at  one  end.  The  curves  marked  as  “Nonlinear”  were  generated  using  a  physically 
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nonlinear  solution.  It  is  obvious  that  the  results  shown  in  Fig.  17,  accounting  for  the 
effect  of  physical  nonlinearity  and  using  the  length  of  the  bimorph  equal  to  2.9  inches  are 
in  excellent  agreement  with  the  experimental  data. 
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Fig.  1.  Piezoelectric  transducers  used  in  ultrasonic  instrumentation  systems  (From: 
http://www.scnscomn.com/lscries.htm  and 

http://www.ndt- 

ed.org/EducationResources/CommunityCollege/Ultrasonics/EquipmentTrans/characterist 

icspt.htm). 
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Fig.  2.  Cylindrical  piezoelectric  rod  and  the  coordinate  system  adopted  in  the  analysis. 


Fig.  3.  The  ratio  of  the  amplitude  obtained  without  accounting  for  physical  nonlinearity 
to  the  physically  nonlinear  counterpart  as  a  function  of  the  nondimensional  frequency  for 
a  PZT-5H  rod.  Damping  is  neglected.  The  electric  field  corresponds  to  O.IMV/m, 
0.3MV/m  and  0.5MV/m  for  cases  1,  2  and  3,  respectively. 
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Fig.  4.  The  ratio  of  the  amplitude  obtained  without  accounting  for  physical  nonlinearity 
to  the  physically  nonlinear  counterpart  as  a  function  of  the  nondimensional  frequency  for 
a  PZN-4.5%PT  rod.  Damping  is  neglected.  The  electric  field  corresponds  to  l.OMV/m, 
1.5MV/m  and  2.0MV/m  for  cases  1,  2  and  3,  respectively. 
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Fig.  5.  The  ratio  of  the  amplitude  of  a  PZT-5H  rod  with  a  physically  linear  material 
behavior  to  the  amplitude  of  the  rod  accounting  for  the  physically  nonlinear  effect  in  the 
presence  of  damping.  The  quality  factor  is  Q  =  65  .  The  half-depth  of  the  rod  is 
h  =  0.1m  .  The  electric  field  corresponds  to  O.IMV/m,  0.3MV/m  and  0.5MV/m  for  cases 
1, 2  and  3,  respectively. 


32 


Fig.  6.  The  ratio  of  the  amplitude  of  a  PZN-4.5%PT  rod  with  a  physically  linear  material 
behavior  to  the  amplitude  of  the  rod  accounting  for  the  physically  nonlinear  effect  in  the 
presence  of  damping.  The  quality  factor  is  Q  =  65  .  The  half-depth  of  the  rod  is 
h  =  0.1  m  .  The  electric  field  corresponds  to  l.OMV/m,  1.5MV/m  and  2.0MV/m  for  cases 
1, 2  and  3,  respectively. 
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Fig.  7.  The  ratio  of  the  amplitude  of  vibrations  of  a  physically  nonlinear  PZT-5H  rod  with 
damping  to  the  amplitude  without  damping  as  a  function  of  the  nondimensional 
frequency.  The  electric  field  is  equal  to  0.5MV/m. 
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Fig.  8.  The  ratio  of  the  amplitude  of  vibrations  of  a  physically  nonlinear  PZN-4.5%PT 
rod  with  damping  to  the  amplitude  without  damping  as  a  function  of  the  nondimensional 
frequency.  The  electric  filed  is  equal  to  2.0MV/m. 
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Fig.  9.  The  ratio  of  the  amplitude  of  vibrations  obtained  using  power  series  to  the 
amplitude  obtained  using  normal  modes  as  a  function  of  the  nondimensional  frequency 
for  a  physically  nonlinear  PZT-5H  rod.  Three  curves  corresponding  to  the  electric  field 
equal  to  O.IMV/m,  0.3MV/m  and  0.5MV/m  practically  coincide. 
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Fig.  10.  The  ratio  of  the  amplitude  of  vibrations  obtained  using  power  series  to  the 
amplitude  obtained  using  normal  modes  as  a  function  of  the  nondimensional  frequency 
for  a  physically  nonlinear  PZN-4.5%PT  rod.  The  electric  field  corresponds  to  l.OMV/m, 
1.5MV/m  and  2.0MV/m  for  cases  1,  2  and  3,  respectively. 
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Fig.  1 1 .  The  ratio  of  the  amplitude  of  motion  corresponding  to  the  second  normal  mode 
to  that  for  the  first  normal  mode  as  a  function  of  the  nondimensional  frequency  for  a 
PZT-5H  rod.  The  curves  for  the  electric  fields  equal  to  O.IMV/m,  0.3MV/m  and 
0.5MV/m  shown  by  solid,  dotted  and  dashed  lines,  respectively,  practically  coincide. 
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Fig.  12.  The  ratio  of  the  amplitude  of  the  motion  corresponding  to  the  second  nonnal 
mode  to  that  for  the  first  normal  mode  as  a  function  of  the  nondimensional  frequency  for 
a  PZN-4.5%PT  rod.  The  electric  field  corresponds  to  l.OMV/m,  1.5MV/m  and  2.0MV/m 
for  cases  1,  2  and  3,  respectively. 
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Fig.  13.  The  nondnnensional  amplitude  of  axial  vibrations  of  a  PZT-5H  rod  (h=5mm)  as 
a  function  of  the  nondimensional  frequency. 


Electric  field:  O.IMV.m;  . 0.3MV/m;  0.5MV/m. 


Fig.  14.  The  nondimensional  amplitude  of  axial  vibrations  of  a  PZT-5H  rod  (h=25mm)  as 
a  function  of  the  nondimensional  frequency. 


Electric  field:  O.IMV.m;  . 0.3MV/m;  0.5MV/m. 
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Fig.  15.  The  nondimensional  amplitude  of  axial  vibrations  of  a  PZN-4.5%PT  rod 
(h=25mm)  as  a  function  of  the  nondimensional  frequency. 

Electric  field:  l.OMV.m;  1.5MV/m;  2.0MV/m. 
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Fig.  16.  Cross-section  of  Model  V  bimorph  (cross  section  along  the  axis  of  the  bimorph). 
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Fig.  17.  Comparison  between  theoretically  predicted  deflections  of  the  tip  of  Model  V 
for  design  2  (linear  and  nonlinear  results)  with  experimental  data  (From  Design  of 
Thunder  Actuators:  Estimate  of  Deformations  and  Analysis  of  a  Discrepancy  between  the 

Analytical  Results  and  Test  Data,  Report  of  Victor  Birman,  LLC  to  QorTek,  Inc.) 
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